In [1]:
#!pip install yellowbrick
#!pip install squarify
#!pip install function_utils
Import des librairies¶
In [89]:
# Import des librairies
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
%matplotlib inline
pd.set_option('display.max_columns', 100)
import missingno as msno
import plotly.graph_objects as go
import plotly.express as px
from sklearn.preprocessing import StandardScaler
from sklearn.preprocessing import MinMaxScaler
from sklearn.metrics import silhouette_score, adjusted_rand_score
from sklearn import metrics
import category_encoders as ce
from plotly.subplots import make_subplots
import datetime
from sklearn.decomposition import PCA
from sklearn.cluster import KMeans
from sklearn.cluster import DBSCAN
from yellowbrick.cluster import KElbowVisualizer
from yellowbrick.cluster import intercluster_distance
from yellowbrick.cluster import SilhouetteVisualizer
import squarify
Import dataset¶
In [3]:
df_olist =pd.read_csv('Dataset/olist_final_cleaned.csv',index_col=0,parse_dates=['order_purchase_timestamp','order_approved_at','order_delivered_carrier_date','order_delivered_customer_date','order_estimated_delivery_date','shipping_limit_date'])
df_olist.info()
In [4]:
df_olist.head()
Out[4]:
In [5]:
df_olist.shape
Out[5]:
Création d'une segmentation RFM¶
- Récence : le nombre de jours depuis la dernière commande
- Fréquence : le nombre de commande passées sur la période étudiée
- Valeur Monétaire : totale du prix des commandes sur la période étudiée
In [6]:
max_date = max(df_olist['order_purchase_timestamp']) + datetime.timedelta(days=1)
In [7]:
rfm_df = df_olist.groupby('customer_unique_id').agg({
'order_purchase_timestamp': lambda x: (max_date - x.max()).days,
'order_id':'count',
'price':'sum'
}).reset_index()
rfm_df.columns = ['customer_unique_id','recency','frequency','monetary']
In [8]:
#snapshot_date = max(df_olist.order_approved_at) + datetime.timedelta(days=1)
In [9]:
#rfm = df_olist.groupby("customer_unique_id").agg({
# "order_approved_at" : lambda x: (snapshot_date - x.max()).days,
# "order_id" : 'count',
# "payment_value" : "mean"
#})
#rfm.columns = ["Recency", "Frequency", "MonetaryValue"]
In [10]:
rfm_df
Out[10]:
In [11]:
rfm_df[rfm_df['frequency'] > 1].count()*100 / rfm_df.shape[0]
Out[11]:
In [12]:
rfm_df.shape
Out[12]:
In [13]:
sns.countplot(y="frequency", data=rfm_df)
plt.title("Répartition de la fréquence d'achat")
plt.show()
La plupart des clients réalise un seul achat
In [14]:
fig, axes = plt.subplots(1, 2, figsize=(15, 5), sharex=True)
fig.suptitle('Description des montants dépensés par les clients')
sns.histplot(ax=axes[0], x="monetary", data=rfm_df, kde=True)
axes[0].set_title("Distribution des montants dépensés")
sns.boxplot(ax=axes[1], x="monetary", data=rfm_df)
axes[1].set_title("Boxplot des montants dépensés")
plt.show()
Les sommes dépensées atteignent un maximum de 14000.
In [15]:
fig, axes = plt.subplots(1, 3, figsize=(15, 5), sharex=True)
fig.suptitle('Distribution des variables')
sns.histplot(ax=axes[0], x="recency", data=rfm_df, kde=True)
axes[0].set_title("Distribution des jours passés depuis la dernière commande")
sns.histplot(ax=axes[1], x="frequency", data=rfm_df, kde=True)
axes[1].set_title("Distribution de la fréquence d'achats des clients")
sns.histplot(ax=axes[2], x="monetary", data=rfm_df, kde=True)
axes[2].set_title("Distribution des sommes dépensées par les clients")
plt.tight_layout()
plt.show()
In [16]:
# Plot RFM distributions
plt.figure(figsize=(15,5))
plt.subplot(1, 3, 1); sns.histplot(rfm_df['recency'])
plt.subplot(1, 3, 2); sns.histplot(rfm_df['frequency'])
plt.subplot(1, 3, 3); sns.histplot(rfm_df['monetary'])
# Show the plot
plt.show()
In [17]:
rfm_df['frequency'].unique()
Out[17]:
Pour la segmentation de notre clientèle à l'aide de K-Means, il sera très important de nous assurer que nous mettons à l'échelle nos données pour centrer la moyenne et les écarts-types. Continuons avec le .qcut() pour notre RFM.
In [18]:
# --Calculate R and F groups--
# Create labels for Recency and Frequency
r_labels = range(4, 0, -1); f_labels = range(1, 5)
# Assign these labels to 4 equal percentile groups
r_groups = pd.qcut(rfm_df['recency'], q=4, labels=r_labels)
# Assign these labels to 4 equal percentile groups
f_groups = pd.qcut(rfm_df['frequency'].rank(method='first'), q=4, labels=f_labels)
# Create new columns R and F
data_process = rfm_df.assign(R = r_groups.values, F = f_groups.values)
data_process.head()
Out[18]:
In [19]:
sns.boxplot(x=data_process['recency'])
Out[19]:
In [20]:
sns.boxplot(x=data_process['frequency'])
Out[20]:
In [21]:
sns.boxplot(x=data_process['monetary'])
Out[21]:
In [22]:
# Create labels for MonetaryValue
m_labels = range(1, 5)
# Assign these labels to three equal percentile groups
m_groups = pd.qcut(data_process['monetary'], q=4, labels=m_labels)
# Create new column M
data_process = data_process.assign(M = m_groups.values)
In [23]:
# Concat RFM quartile values to create RFM Segments
def join_rfm(x): return str(x['R']) + str(x['F']) + str(x['M'])
data_process['RFM_Segment_Concat'] = data_process.apply(join_rfm, axis=1)
rfm = data_process
rfm.head()
Out[23]:
In [24]:
# Count num of unique segments
rfm_count_unique = rfm.groupby('RFM_Segment_Concat')['RFM_Segment_Concat'].nunique()
print(rfm_count_unique.sum())
In [25]:
# Calculate RFM_Score
rfm['RFM_Score'] = rfm[['R','F','M']].sum(axis=1)
print(rfm['RFM_Score'].head())
In [26]:
rfm
Out[26]:
Analyzing RFM Segmentation Let’s delve into few interesting segments:
Champions are your best customers, who bought most recently, most often, and are heavy spenders. Reward these customers. They can become early adopters for new products and will help promote your brand.
Potential Loyalists are your recent customers with average frequency and who spent a good amount. Offer membership or loyalty programs or recommend related products to upsell them and help them become your Loyalists or Champions.
New Customers are your customers who have a high overall RFM score but are not frequent shoppers. Start building relationships with these customers by providing onboarding support and special offers to increase their visits.
At Risk Customers are your customers who purchased often and spent big amounts, but haven’t purchased recently. Send them personalized reactivation campaigns to reconnect, and offer renewals and helpful products to encourage another purchase.
Can’t Lose Them are customers who used to visit and purchase quite often, but haven’t been visiting recently. Bring them back with relevant promotions, and run surveys to find out what went wrong and avoid losing them to a competitor.
In [27]:
# Define rfm_level function
def rfm_level(df):
if df['RFM_Score'] >= 9:
return 'Can\'t Lose Them'
elif ((df['RFM_Score'] >= 8) and (df['RFM_Score'] < 9)):
return 'Champions'
elif ((df['RFM_Score'] >= 7) and (df['RFM_Score'] < 8)):
return 'Loyal'
elif ((df['RFM_Score'] >= 6) and (df['RFM_Score'] < 7)):
return 'Potential'
elif ((df['RFM_Score'] >= 5) and (df['RFM_Score'] < 6)):
return 'Promising'
elif ((df['RFM_Score'] >= 4) and (df['RFM_Score'] < 5)):
return 'Needs Attention'
else:
return 'Require Activation'
# Create a new variable RFM_Level
rfm['RFM_Level'] = rfm.apply(rfm_level, axis=1)
# Print the header with top 5 rows to the console
rfm.head()
Out[27]:
In [29]:
# Calculate average values for each RFM_Level, and return a size of each segment
rfm_level_agg = rfm.groupby('RFM_Level').agg({
'recency': 'mean',
'frequency': 'mean',
'monetary': ['mean', 'count']
}).round(1)
# Print the aggregated dataset
print(rfm_level_agg)
À partir de là, nous pouvons voir qu'un grand pourcentage de nos clients se situent dans les niveaux supérieurs de RFM. Le magasin doit faire quelque chose de bien pour maintenir sa fidélité !
Potential - potentiel élevé pour entrer dans nos segments de clients fidèles
Promising — montre des signes prometteurs avec la quantité et la valeur de leur achat, mais cela fait un moment qu'ils n'ont pas acheté.
Needs Attention — a fait un achat initial mais ne l'a pas vu depuis.
Require Activation — Les moins performants de notre modèle RFM.
In [30]:
rfm_level_agg
Out[30]:
In [31]:
rfm_level_agg.columns = rfm_level_agg.columns.droplevel()
rfm_level_agg.columns = ['RecencyMean','FrequencyMean','MonetaryMean', 'Count']
In [32]:
rfm_level_agg = rfm_level_agg.reset_index()
In [33]:
#Create our plot and resize it.
fig = plt.gcf()
ax = fig.add_subplot()
#fig.set_size_inches(16, 9)
squarify.plot(sizes=rfm_level_agg['Count'],
label=rfm_level_agg['RFM_Level'], color=["red","cyan","blue", "grey", "navy","yellow","green"], alpha=.6 )
plt.title("RFM Segments",fontsize=18,fontweight="bold")
plt.show()
In [34]:
fig = px.scatter_3d(rfm, x='recency',
y='frequency', z='monetary',
color='RFM_Level')
fig.show()
In [35]:
pie = rfm.RFM_Level.value_counts()
plt.figure(figsize=(12, 8))
plt.pie(pie, labels= pie.index, autopct='%1.1f%%')
plt.title('repartion des classes RFM')
plt.show()
On va concatener les 2 datasets
In [36]:
df_merged = df_olist.merge(rfm, how="left", on="customer_unique_id")
In [37]:
df_merged.boxplot('payment_value', by='RFM_Level', figsize=(12, 8))
plt.xticks(rotation=90)
plt.title('prix commande par classe RFM')
plt.show()
In [38]:
df_merged.head()
Out[38]:
In [39]:
fig = px.scatter(df_merged, x=df_merged.recency,
y=df_merged.frequency,
color=df_merged.customer_state,
size=df_merged.monetary,
title="Distribution of Frequency and Recency of purchase as per Customer State<br><sup>Size of Circles represent the purchase monetary</sup>")
fig.show()
In [40]:
fig = px.scatter(df_merged, x=df_merged.monetary,
y=df_merged.frequency,
color=df_merged.customer_state,
size=df_merged.recency,
title="Distribution of Frequency and Moneytary purchase as per Customer State<br><sup>Size of Circles represent the recency of purchase</sup>")
fig.show()
On va enlever les outliers de par rapport au monetary
In [41]:
#On supprime les valeurs aberrantes
def remove_outlier_iqr(df_in, column):
'''function that gives IQR,lower and upper whisker'''
#On calcule Q1 et Q3
q1, q3 = np.quantile(df_in[column], 0.25), np.quantile(df_in[column], 0.95)
#On calcule l'écart interquartile (IQR)
IQR = q3 - q1
#On calcule la borne inférieure à l'aide du Q1 et de l'écart interquartile
borne_inf = q1-1.5*IQR
#On calcule la borne supérieure à l'aide du Q3 et de l'écart interquartile
borne_sup = q3 +1.5*IQR
#On calcul les valeurs à l'intérieur de la borne inférieure et supérieure
df1 = df_in[df_in[column] > borne_sup]
df2 = df_in[df_in[column] < borne_inf]
''' Removing the Outliers '''
df_out = df_in[~((df_in[column] < borne_inf) |(df_in[column] > borne_sup))]
print('Total number of outliers are ', df1.shape[0]+ df2.shape[0])
return df_out
In [42]:
df_stats = remove_outlier_iqr(df_merged, "monetary")
In [43]:
df_stats = remove_outlier_iqr(df_merged, "frequency")
In [44]:
df_stats.shape
Out[44]:
In [45]:
sns.boxplot(x=df_stats['monetary'])
Out[45]:
In [46]:
sns.boxplot(x=df_stats['frequency'])
Out[46]:
In [47]:
df_stats.head()
Out[47]:
In [48]:
df_stats.to_csv("Dataset/olist_df_stats.csv")
In [49]:
recency = np.log(df_stats['recency'])
frequency = np.log(df_stats["frequency"])
monetary = np.log(df_stats["monetary"])
In [50]:
df = df_stats.copy()
In [51]:
df.loc[:,"recency"] = recency
df.loc[:,"frequency"] = frequency
df.loc[:,"monetary"] = monetary
df.reset_index(inplace=True)
In [52]:
fig, axes = plt.subplots(1, 3, figsize=(15, 6), sharex=True)
fig.suptitle('Distribution des variables après transformation logarithme')
sns.histplot(ax=axes[0], x="recency", data=df, kde=True)
axes[0].set_title("Distribution des jours passés depuis la dernière commande")
sns.histplot(ax=axes[1], x="frequency", data=df, kde=True)
axes[1].set_title("Distribution de la fréquence d'achats des clients")
sns.histplot(ax=axes[2], x="monetary", data=df, kde=True)
axes[2].set_title("Distribution des sommes dépensées par les clients")
plt.tight_layout()
plt.show()
On va ajouter d'autres variables pour la segmentation. On va ajouter le nombre de jours entre la commande et la réception de cette dernière, l'avis du client et la catégorie de produits commandée.In [53]:
df["nb_days_to_deliver"] = df["order_delivered_customer_date"] - df["order_approved_at"]
df["nb_days_estimated_to_deliver"] = df["order_estimated_delivery_date"] - df["order_approved_at"]
add_cols = ["nb_days_to_deliver","nb_days_estimated_to_deliver"]
for cols in add_cols:
df[cols] = df[cols].dt.days
df.drop(["order_approved_at",
"order_delivered_customer_date",
"order_estimated_delivery_date"], axis=1, inplace=True)
In [54]:
df[df["nb_days_to_deliver"] < 0]
Out[54]:
In [55]:
df = df.loc[df["nb_days_to_deliver"] > 0]
In [56]:
df[df["nb_days_estimated_to_deliver"] < 0]
Out[56]:
In [57]:
df = df.loc[df["nb_days_estimated_to_deliver"] > 0]
In [58]:
df.head()
Out[58]:
In [59]:
df.to_csv("Dataset/data_stability.csv")
In [60]:
df.boxplot('nb_days_to_deliver', by='RFM_Level', figsize=(12, 8))
plt.xticks(rotation=90)
plt.title('Nombre de jour de livraison commande par classe RFM')
plt.show()
In [61]:
df_segment = df[["customer_unique_id","recency","frequency","monetary"]]
In [62]:
df_segment.set_index("customer_unique_id", inplace=True)
In [63]:
df_segment.head()
Out[63]:
Test avec uniquement les variables RFM¶
Standardisation des données
In [64]:
scaler = StandardScaler()
In [65]:
scaler.fit(df_segment)
rfm_normalized = scaler.transform(df_segment)
Réduction de dimensions¶
La réduction de dimension a pour but de faciliter la visualisation de nos données en les ramenant dans un espace en 2 dimensions.
In [66]:
# PCA Pipeline
pca = PCA(svd_solver='full')
pca.fit(rfm_normalized)
data_pca = pca.transform(rfm_normalized)
In [67]:
# Explained variance
varexpl = pca.explained_variance_ratio_*100
# Plot of cumulated variance
plt.figure(figsize=(12,8))
plt.bar(np.arange(len(varexpl))+1, varexpl)
cumSumVar = varexpl.cumsum()
plt.plot(np.arange(len(varexpl))+1, cumSumVar,c="red",marker='o')
valid_idx = np.where(cumSumVar >= 95)[0]
min_plans = valid_idx[cumSumVar[valid_idx].argmin()]+1
plt.xlabel("rang de l'axe d'inertie")
plt.xticks(np.arange(len(varexpl))+1)
plt.ylabel("pourcentage d'inertie")
plt.title("{}% de la variance totale est expliquée"\
" par les {} premiers axes".format(95,
min_plans))
plt.show(block=False)
On conserve les 3 axes principaux pour expliquer la variance à 95%.
In [68]:
def display_circles(pcs, n_comp, pca, axis_ranks, labels=None, label_rotation=0, lims=None):
for d1, d2 in axis_ranks: # On affiche les 3 premiers plans factoriels, donc les 6 premières composantes
if d2 < n_comp:
# initialisation de la figure
fig, ax = plt.subplots(figsize=(10,10))
# détermination des limites du graphique
if lims is not None :
xmin, xmax, ymin, ymax = lims
elif pcs.shape[1] < 30 :
xmin, xmax, ymin, ymax = -1, 1, -1, 1
else :
xmin, xmax, ymin, ymax = min(pcs[d1,:]), max(pcs[d1,:]), min(pcs[d2,:]), max(pcs[d2,:])
# affichage des flèches
# s'il y a plus de 30 flèches, on n'affiche pas le triangle à leur extrémité
if pcs.shape[1] < 30 :
plt.quiver(np.zeros(pcs.shape[1]), np.zeros(pcs.shape[1]),
pcs[d1,:], pcs[d2,:],
angles='xy', scale_units='xy', scale=1, color="grey")
# (voir la doc : https://matplotlib.org/api/_as_gen/matplotlib.pyplot.quiver.html)
else:
lines = [[[0,0],[x,y]] for x,y in pcs[[d1,d2]].T]
ax.add_collection(LineCollection(lines, axes=ax, alpha=.1, color='black'))
# affichage des noms des variables
if labels is not None:
for i,(x, y) in enumerate(pcs[[d1,d2]].T):
if x >= xmin and x <= xmax and y >= ymin and y <= ymax :
plt.text(x, y, labels[i], fontsize='14', ha='center', va='center', rotation=label_rotation, color="blue", alpha=0.5)
# affichage du cercle
an = np.linspace(0, 2 * np.pi, 100) # Add a unit circle for scale
plt.plot(np.cos(an), np.sin(an))
plt.axis('equal')
# définition des limites du graphique
plt.xlim(xmin, xmax)
plt.ylim(ymin, ymax)
# affichage des lignes horizontales et verticales
plt.plot([-1, 1], [0, 0], color='grey', ls='--')
plt.plot([0, 0], [-1, 1], color='grey', ls='--')
# nom des axes, avec le pourcentage d'inertie expliqué
plt.xlabel('F{} ({}%)'.format(d1+1, round(100*pca.explained_variance_ratio_[d1],1)))
plt.ylabel('F{} ({}%)'.format(d2+1, round(100*pca.explained_variance_ratio_[d2],1)))
plt.title("Cercle des corrélations (F{} et F{})".format(d1+1, d2+1))
plt.show(block=False)
In [70]:
# Principal component space
pcs = pca.components_
display_circles(pcs, 3, pca, [(0, 1)], labels=np.array(df_segment.columns))
K-Means Clustering¶
K-Means se base sur des calculs de distance entre les points de notre jeu de données et un point nommé centroïde.
Le KElbowVisualizer implémente la méthode de "coude" pour sélectionner le nombre optimal de clusters en ajustant le modèle avec une plage de valeurs pour 𝐾. Si le graphique linéaire ressemble à un bras, alors le « coude » (le point d'inflexion sur la courbe) est une bonne indication que le modèle sous-jacent s'adapte le mieux à ce point. Dans le visualiseur, "coude" sera annoté par une ligne en pointillés.
Le KElbowVisualizer affiche également le temps de calcul de l'entrainement du modèle de clustering par 𝐾 sous la forme d'une ligne verte en pointillés
In [71]:
# Instantiate the clustering model and visualizer
model = KMeans()
visualizer = KElbowVisualizer(
model, k=(2,10), metric='calinski_harabasz', timings=True
)
visualizer.fit(data_pca) # Fit the data to the visualizer
visualizer.show() # Finalize and render the figure
Out[71]:
In [72]:
# Instantiate the clustering model and visualizer
model = KMeans()
visualizer = KElbowVisualizer(model, k=(2,12))
visualizer.fit(data_pca) # Fit the data to the visualizer
visualizer.show() # Finalize and render the figure
Out[72]:
Selon la méthode du coude basée sur le score de distortion (somme moyenne des carrés des distances aux centres), le nombre de cluster idéal est de 4.
On va donc entraîner notre modèle avec ce paramètre.
In [73]:
# Fitting KMeans
km = KMeans(n_clusters=4, max_iter=50)
In [74]:
# Prediction
label = km.fit_predict(data_pca)
centroids = km.cluster_centers_
u_labels = np.unique(label)
# Graphical representation
plt.figure(figsize=(10, 10))
for i in u_labels:
plt.scatter(data_pca[label==i, 0], data_pca[label==i, 1], label = i)
plt.scatter(centroids[:, 0], centroids[:,1], s=80, alpha=0.8, color='k')
plt.legend()
plt.title('Représentation des clusters du K-Means')
plt.show()
Silhouette Visualizer¶
L'analyse de silhouette est utilisée pour évaluer la densité et la séparation entre les clusters. Le score est calculé en faisant la moyenne du coefficient de silhouette pour chaque échantillon, qui est calculé comme la différence entre la distance moyenne des clusters et la distance moyenne pour chaque échantillon, normalisée par la valeur maximale. Cela produit un score compris entre -1 et +1, où les scores proches de +1 indiquent une séparation élevée et les scores proches de -1 indiquent que les échantillons ont peut-être été affectés au mauvais groupe.
Pour vérifier si le Kmeans clustering fonctionne, nous allons utiliser SilhouetteVisualizer pour afficher le coefficient de silhouette pour un échantillonage de chaque cluster. Cela permet de visualiser la densité et la séparation des clusters.
In [76]:
fig, ax = plt.subplots(2, 2, figsize=(15,8))
for i in [2, 3, 4, 5]:
'''
Create KMeans instance for different number of clusters
'''
k_means_model = KMeans(n_clusters=i, init='k-means++', n_init=10, max_iter=100, random_state=42)
q, mod = divmod(i, 2)
'''
Create SilhouetteVisualizer instance with KMeans instance
Fit the visualizer
'''
visualizer = SilhouetteVisualizer(k_means_model, colors='yellowbrick', ax=ax[q-1][mod])
visualizer.fit(df_segment)
Le graphique ci-dessus est une comparaison entre 4 scénarios de clustering différents, où chaque scénario représente un certain nombre de clusters, 2 à 5 respectivement. Au fur et à mesure que le nombre de clusters augmente, nous voyons ce qui suit :
le score Silhouette moyen de tous les clusters diminue, représenté par la ligne pointillée verticale, ce qui signifie que nous n'optimisons pas à mesure que nous augmentons le nombre de clusters le saut du score moyen est clair de 2 clusters à 3 clusters (jusqu'à ~ 0,16 de ~ 0,35) indiquant que 2 est la valeur optimale
à mesure que nous augmentons les clusters, plus de points ont tendance à avoir un score Silhouette négatif confirmant que les clusters commencent à devenir moins homogènes on remarque aussi que pour k = 2, le premier cluster est presque trois fois plus grand que la taille (nombre d'instances) du second cluster
In [77]:
from yellowbrick.cluster import SilhouetteVisualizer
visualizer = SilhouetteVisualizer(km, colors='yellowbrick')
visualizer.fit(data_pca)
visualizer.show()
Out[77]:
In [78]:
# Kmeans labels
df_segment["cluster"] = km.labels_
df_segment
Out[78]:
Analyse des clusters¶
In [79]:
fig = px.scatter_3d(df_segment, x=df_segment["recency"], y=df_segment["frequency"], z=df_segment["monetary"], color=df_segment["cluster"])
fig.show()
In [80]:
sns.boxplot(x='cluster',y='monetary', data=df_segment)
Out[80]:
In [81]:
sns.boxplot(x='cluster',y='frequency', data=df_segment)
Out[81]:
In [82]:
sns.boxplot(x='cluster',y='recency', data=df_segment)
Out[82]:
In [83]:
# Cluster 0
In [84]:
cluster0 = df_segment[df_segment['cluster'].isin([0])]
In [85]:
customers_0 = cluster0.index
df0 = df.loc[df['customer_unique_id'].isin(customers_0)]
df0.shape
Out[85]:
In [86]:
df0.head()
Out[86]:
In [87]:
product_cluster0 = df0.groupby('product_category_name').agg({'payment_value': lambda payment_value: payment_value.sum()})
DBSCAN¶
In [90]:
# Compute DBSCAN
db = DBSCAN(eps = 1, min_samples=5,n_jobs=-1).fit(data_pca)
core_samples_mask = np.zeros_like(db.labels_, dtype=bool)
core_samples_mask[db.core_sample_indices_] = True
labels2D = db.labels_
# Number of clusters in labels, ignoring noise if present.
n_clusters_ = len(set(labels2D)) - (1 if -1 in labels2D else 0)
n_noise_ = list(labels2D).count(-1)
print('2D Estimated number of clusters: %d' % n_clusters_)
print('2D Estimated number of noise points: %d' % n_noise_)
print("2D Silhouette Coefficient: %0.3f" % metrics.silhouette_score(data_pca, labels2D))
In [91]:
db1 = DBSCAN(eps=0.4, min_samples=3).fit(data_pca)
core_samples_mask = np.zeros_like(db1.labels_, dtype=bool)
core_samples_mask[db1.core_sample_indices_] = True
labels2D2 = db1.labels_
# Number of clusters in labels, ignoring noise if present.
n_clusters_ = len(set(labels2D2)) - (1 if -1 in labels2D2 else 0)
n_noise_ = list(labels2D2).count(-1)
print('2D Estimated number of clusters: %d' % n_clusters_)
print('2D Estimated number of noise points: %d' % n_noise_)
print("2D Silhouette Coefficient: %0.3f" % metrics.silhouette_score(data_pca, labels2D2))
In [92]:
db2 = DBSCAN(eps=0.1, min_samples=10).fit(data_pca)
core_samples_mask = np.zeros_like(db2.labels_, dtype=bool)
core_samples_mask[db2.core_sample_indices_] = True
labels2D3 = db2.labels_
# Number of clusters in labels, ignoring noise if present.
n_clusters_ = len(set(labels2D3)) - (1 if -1 in labels2D3 else 0)
n_noise_ = list(labels2D3).count(-1)
print('2D Estimated number of clusters: %d' % n_clusters_)
print('2D Estimated number of noise points: %d' % n_noise_)
print("2D Silhouette Coefficient: %0.3f" % metrics.silhouette_score(data_pca, labels2D3))
In [93]:
# Création d'un dataset avec features et labels
df_dbscan = df.copy()
df_dbscan['labels'] = db.labels_
In [94]:
stat = pd.pivot_table(df_dbscan,
values=df_dbscan.columns,
index=['labels', 'customer_unique_id'])
for x in df_dbscan['labels'].value_counts().index:
display(stat.loc[x, :].describe())
In [95]:
df_dbscan.columns
Out[95]:
In [96]:
df_dbsc = df_dbscan[['payment_type', 'payment_value', 'review_score', 'price','product_category_name',
'recency', 'frequency', 'monetary', 'RFM_Level', 'nb_days_to_deliver',
'nb_days_estimated_to_deliver', 'labels']]
In [97]:
# Représentation graphique des features pour chaque cluster
for col in df_dbsc.columns.to_list():
plt.figure(figsize=(15, 7))
sns.boxplot(y=col, x='labels', data=df_dbsc,
showfliers=False)
plt.show()
In [98]:
# Analyse des différentes catégories dans les labels
index_tot = [df_dbsc[df_dbsc['labels'] == x].index
for x in df_dbsc['labels'].value_counts().index]
plt.figure(figsize=(20, 12))
for x in range(len(index_tot)):
order = df_dbsc.loc[index_tot[x], 'payment_type'].value_counts()
order_hue = order.index
plt.subplot(2, len(index_tot)/2, x+1)
sns.countplot(y=df_dbsc.loc[index_tot[x], 'payment_type'],
order=order_hue,
palette='Blues_r')
plt.title(f"Cluster {x}", fontsize=20)
plt.figure(figsize=(20, 12))
for x in range(len(index_tot)):
order = df_dbsc.loc[index_tot[x], 'product_category_name'].value_counts()
order_hue = order.index
plt.subplot(2, len(index_tot)/2, x+1)
sns.countplot(y=df_dbsc.loc[index_tot[x], 'product_category_name'],
order=order_hue,
palette='Blues_r')
plt.title(f"Cluster {x}", fontsize=20)
In [99]:
def display_factorial_planes(X_projected, n_comp, pca, axis_ranks, hue_serie=None,
xlimit=None, ylimit=None, labels=None,
alpha=1, illustrative_var=None):
"""Tracés des projections des individus dans les plans factoriels"""
import matplotlib
matplotlib.rcdefaults()
for d1,d2 in axis_ranks:
if d2 < n_comp:
# initialisation de la figure
fig = plt.figure(figsize=(7, 6))
# affichage des points
if illustrative_var is None:
sns.scatterplot(x=X_projected[:, d1],
y=X_projected[:, d2],
hue=hue_serie,
alpha=alpha,
palette='tab10')
else:
illustrative_var = np.array(illustrative_var)
for value in np.unique(illustrative_var):
selected = np.where(illustrative_var == value)
sns.scatterplot(x=X_projected[selected, d1],
y=X_projected[selected, d2],
hue=hue_serie,
alpha=alpha,
label=value,
palette='tab10')
plt.legend()
# affichage des labels des points
if labels is not None:
for i,(x,y) in enumerate(X_projected[:, [d1, d2]]):
plt.text(x, y, labels[i],
fontsize='14', ha='center',va='center')
# détermination des limites du graphique
boundary = np.max(np.abs(X_projected[:, [d1, d2]])) * 1.1
if (xlimit, ylimit) == (None, None) :
plt.xlim([-boundary, boundary])
plt.ylim([-boundary, boundary])
else:
plt.xlim([-xlimit, xlimit])
plt.ylim([-ylimit, ylimit])
# affichage des lignes horizontales et verticales
plt.plot([-100, 100], [0, 0], color='grey', ls='--')
plt.plot([0, 0], [-100, 100], color='grey', ls='--')
# nom des axes, avec le pourcentage d'inertie expliqué
plt.xlabel('F{} ({}%)'.format(d1+1, round(100*pca.explained_variance_ratio_[d1], 1)))
plt.ylabel('F{} ({}%)'.format(d2+1, round(100*pca.explained_variance_ratio_[d2], 1)))
plt.title("Projection des individus (sur F{} et F{})".format(d1+1, d2+1))
plt.savefig(f"{d1}_{d2}_factorial_plane.png")
plt.show(block=False)
In [100]:
# Projection des individus
display_factorial_planes(data_pca, 3, pca,
[(0, 1), (1, 2)], hue_serie=df_dbsc['labels'])
plt.show()
Stabilité des 2 clustering¶
In [101]:
from sklearn.base import clone
from sklearn.utils import check_random_state
rng = np.random.RandomState(1)
def cluster_stability(X, est, n_iter=20, random_state=None):
labels = []
indices = []
for i in range(n_iter):
# draw bootstrap samples, store indices
sample_indices = rng.randint(0, X.shape[0], X.shape[0])
indices.append(sample_indices)
est = clone(est)
if hasattr(est, "random_state"):
# randomize estimator if possible
est.random_state = rng.randint(1e5)
X_bootstrap = X[sample_indices]
est.fit(X_bootstrap)
# store clustering outcome using original indices
relabel = -np.ones(X.shape[0], dtype=np.int)
relabel[sample_indices] = est.labels_
labels.append(relabel)
scores = []
for l, i in zip(labels, indices):
for k, j in zip(labels, indices):
# we also compute the diagonal which is a bit silly
in_both = np.intersect1d(i, j)
scores.append(adjusted_rand_score(l[in_both], k[in_both]))
return np.mean(scores)
In [102]:
km_stability = []
cluster_range = range(2, 10, 1)
for n_clusters in cluster_range:
print(n_clusters)
km = KMeans(n_clusters=n_clusters, n_init=10, init="random")
km_stability.append(cluster_stability(data_pca, km))
In [103]:
db_stability = []
n_clusters_db = []
for eps in np.linspace(.2, 1, 10):
print(eps)
db_stability.append(cluster_stability(data_pca, DBSCAN(eps=eps)))
n_clusters_db.append(len(np.unique(DBSCAN(eps=eps).fit(data_pca).labels_)))
In [104]:
cluster_range = range(2, 10, 1)
plt.plot(cluster_range, km_stability, label="k-means")
plt.plot(n_clusters_db, db_stability, label="DBSCAN")
for eps, n_clusters, stability in zip(np.linspace(.2, 1, 15), n_clusters_db, db_stability):
plt.text(n_clusters, stability, "{:.2f}".format(eps))
plt.legend()
plt.xlabel("n_clusters")
plt.ylabel("stability")
Out[104]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]: